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Abstract 

Classical vector analysis is the predominant formalism used by 

engineers of computational electromagnetism, despite the fact that 
manifold as a theoretical concept has existed for a century. This pa- 
per discusses the benefits of manifolds over the traditional approach 
in practical problems of modelling. With a structural approach, it 
outlines the role and interdependence of coordinate systems, met- 
ric, constitutive equations, and fields, and relates them to practical 
problems of quasi-static computational electromagnetics: mesh gener- 
ation, open-boundary problems, and electromagnetic-mechanical cou- 
pled problems involving motion and deformation. The proposed pro- 
cedures also imply improvements to the flexibility of the modelling 
software. 

1 Introduction 

The objective of solver software systems for electromagnetic boundary value 
problems (BVP) is to solve an arbitrary problem from a specific problem 
class. However, this is not possible in practice: One of the restricting factors 
is the difficulty of the mesh generation. This is related to the fact that 
to pose a BVP, the user needs to define a geometry of the BVP domain. 
The geometry is described with a coordinate system and certain choices of 
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coordinate systems may cause too much round-off error with ffoating point 
arithmetic. 

Often the problem geometries lend themselves to a particular coordinate 
system, and therefore the software systems usually offer a list of coordinate 
systems to choose from. To offer a user interface with a list of coordinate 
systems, many structures in the software are fixed: For example, the soft- 
ware often has one intrinsic coordinate system, usually Cartesian, which it 
uses in all calculations. Furthermore, the basic material parameters, such as 
the permittivity of a vacuum, as well as the metric (inner product) in the 
intrinsic coordinate system are pre-decided. When the user gives the model's 
geometry, for example, in the spherical coordinates, the software maps them 
to the intrinsic Cartesian system. 

Frequently, the floating point problems leave the modeler without any 
mesh at all and prevent the solution of the BVP: Practical computations 
with the intrinsic coordinate system requires that the coordinate numbers 
are represented with floating point numbers. This makes mesh generation 
particularly challenging for BVPs with details of the domain varying signifl- 
cantly in size. For example, consider the power line on top of flgure [Tj The 
thickness of the cables and the structures in the supports are much smaller in 
scale than the distance between the two supports. Meshes are difficult to gen- 
erate because the floating point representation gives us labels to only a flnite, 
although large, number of points. For example, as a floating point number 
grows in magnitude, the distance between points separable by numeric labels 
increases, and, at some point, round-off errors become big enough to make 
labels of distinct points coincide unintentionally. Generally, mesh generation 
is based on geometric predicates [Ij, some of which are difficult to calculate 
with floating point numbers. Mesh generation becomes hard or even impos- 
sible with double precision floating point numbers when the scale variation, 
the ratio of the largest and the smallest dimensions, is in the order of 10^ or 
higher. This estimate for scale variation is based on the experiences with the 
BVP presented in the section 7. 

The mesh generation problems arise from the way the BVP was described 
in a particular coordinate system. Is it possible to alleviate these problems? 
With the help of modern mathematics, manifolds and differential geometry, 
the role of mathematical structures on the description of a BVP can be out- 
lined. In this paper we strongly promote a structural approach. Particularly, 
we will show that for a given BVP, the elements of the triplet {coordinate 
system, metric, material parameters} are interdependent, and the triplet con- 
stitutes an entity. In fact, there is a whole class of these triplets for each BVP: 
Each triplet describes a geometry of the BVP domain and the constitutive 
equations, which are needed to deflne an electromagnetic BVP. Speciflcally, 
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Figure 1: Power line in two coordinate system. 



the numerical values of material parameters are respect to given metric. For 
example, the numerical value of the permittivity of vacuum is different with 
respect to meters than with respect to inches. Although the metric appears 
in the triplets, the electromagnetic theory essentially concerns the works re- 
lated to displacements of particles and the total energy of systems. The role 
of metric is central only in the classical view — vector analysis — whereas the 
modern view of differential geometry reveals the independence of the elec- 
tromagnetic theory from a particular metric. This has been known for long 
[2] [3] [4] [5] [6] , but not that well recognized in engineering. Hence, the modeller 
may select the metric at will. This is much greater a freedom in the context 
of manifold than a redefinition of the inner product in the single codomain 
vector space of a vector field in the classical view. We use the modern view 
to show how to formulate electromagnetic problems with the classical vector 
formalism in different metrics. 

From the software design point of view, there is no need to fix the whole 
triplet in a software system. Just like the software system is not intended 
to solve a particular problem, but any problem from a class of problems, 
the solver software should provide the possibility to use any triplet from the 
triplet class. It should not be even restricted to a single triplet, because 
manifold promotes local coordinate systems that provide further flexibility. 
This will help overcome floating point problems in the calculations of the 
geometric predicates and thus will enlarge the class of practically solvable 
problems: Our approach suggests giving up the metric and the coordinate 
system strongly suggested by human intuition. For example, if the mesh 
generation fails with the coordinate system in the top of figure [T| one could 
use the coordinate system in the bottom of the same figure, when the calcu- 
lations of the predicates are easier. This time it is more likely that one gets 
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a mesh and thus a numerical solution. However, there is usually a price to 
pay: numerical error can be large because the mesh quality may be poor. 

The problems of mesh generation are just one example of practical chal- 
lenges in numerics that can be engaged with tools provided by manifold. The 
problems with motion and deformations, such as the deformation effects of 
magnetostriction, can be modelled with a single mesh: motions and defor- 
mations can be modelled as changes in the metric. This can speed up the 
solution process substantially. Finally, the solutions for open boundary prob- 
lems, that are solved by finite domains with so-called infinite elements in the 
boundary region, are an example of an application where a convenient triplet 
is chosen for solution. 

2 Structural Background 

We need many mathematical structures to pose and describe an electromag- 
netic EVP. The starting point of this description is a point set that is the 
EVP domain. The points of this set correspond to observations on distance 
measurements. Then we will add one structure at a time to this point set, 
such that finally we have all that is needed for the description. 

2.1 Metric Space 

At first we need to specify a point set for the EVP domain. The points in the 
domain correspond to the points observed in reality, in the following sense: 
The observed points are considered distinct if and only if a nonzero distance 
can be measured between them. We denote this abstract point set by M. Ey 
the way the set M was constructed, it is naturally endowed with the structure 
of metric space [7]. The distance between a given pair of points refers to an 
observation of lengths in multiples of some reference object called the rigid 
body. (Observe here that numerical values of distances depend on the chosen 
rigid body but distances can be discussed without a specified reference body.) 
Furthermore, the metric induces the structure of topological space for M [7]. 
The dimension of M, which is two or three, must be the same at every point 
of space M. 

2.2 Manifolds 

For computational analysis, we need to characterize displacements arithmeti- 
cally and thus must parameterize the point set. In the classical view this is 
done with coordinate systems, which assign unique coordinates — that is, a 
tuple of real numbers — to each point of M . Formally this means a mapping 
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of type M M". Furthermore, M" is endowed with the norm topology and 
the mapping is homeomorphic, i.e. continuous bijection with continuous in- 
verse. However, the abstract point set M is often ignored completely and 
the points of M are identified with the coordinate numbers. Furthermore, 
often a metric connotation is given to the coordinate numbers. Thus the 
domain of the EVP is taken to be some particular subset of M". Yet the co- 
ordinates are just labels for the points of M, such that they allow arithmetic 
and need not carry metric information. The coordinates can be changed with 
diffeomorphic mappings of type M" . 

Instead of a particular coordinate system, the modern view uses manifolds 
as domains for EVPs. In manifolds, the abstract point set M is the primary 
object, which is not identified with any particular set of coordinate numbers. 
The coordinate systems, or charts, are considered secondary, which is natural 
because we can choose the charts in many ways, and the real numbers as- 
signed to each point are not canonical. In particular, the physics — work and 
energy — does not depend at all on the chosen chart. Furthermore, M is not 
expected to be coverable by a single chart, as was the case for coordinate sys- 
tems: Charts are homeomorphic mappings from the open connected subsets 
of M to the open subsets of M". Moreover, every point of M is required to be 
in the domain of some chart. Any topological space, that can be covered with 
charts, is called a topological manifold [8J. Thus the structure of topological 
manifold emphasizes the existence of charts, rather than a particular chart. 

Topological manifold is not enough for physics, because the differentia- 
bility of functions needs an extra structure. Let us require a collection of 
charts, that are bound together by their transition maps, or change-of-chart 
maps: The differentiability of a real-valued function h defined on M is de- 
fined indirectly through the differentiability of the functions h o /~^, where 
/ is a chart. The chart can be changed with the transition map f o as 
follows: {h o /^^) o (/ o g^^) = ho g^^. We do not want the n-times differen- 
tiability of h to depend on the choice of chart, and therefore require that the 
transition map be (at least) n-times different iable. If this holds, the charts 
are considered equivalent. This relation constitutes an equivalence class, a 
difjerentiable structure, of charts. The set M together with a differentiable 
structure is a differentiable manifold or here just manifold for short [8][9]|10]. 

2.3 Standard Parameterization 

The rigid body can be used as a pair of dividers to specify spheres in M. 
If the image of every sphere in M under a chart is a sphere also in 
when measured with the Euclidean 2-norm, we call the chart a standard 
parameterization. For example, the chart / in figure |2] is a standard pa- 
rameterization. All standard parameterizations are identical up to scaling, 
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Figure 2: Standard parameterization: The line drawing on top refers to a real 
coaxial cable. The rigid body we chose gave us a pair of dividers and enabled 
us to specify spheres. Chart / from M to is a standard parameterization, 
for the images of the spheres in M are also spheres in M"^ in the sense of 
Euclidean 2-norm. Consequently, chart g is not a standard parameterization. 




Figure 3: Tangent space and tangent vector. A tangent vector v of the 
tangent space of the point p in manifold M. 
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rotation, and translation (plus reflection). The standard parameterizations 
are important because in literature the numerical values of the material pa- 
rameters are given in a specific class of standard parameterizations. This is 
natural because the numerical values of the material parameters are mean- 
ingful only with respect to a given metric, and in practice, the given metric 
is always based on a length unit system in which some rigid body is chosen 
for reference. 

2.4 Tangent Spaces 

Charts enable us to talk about displacements in small neighborhoods, i.e. 
virtual translations: the possibility to parameterize M implies a local n- 
dimensional vector space at each point of M. The vector space is called the 
tangent space, and the elements of this space are tangent vectors [8][9]il0j[llJ. 
Virtual translations can be recognized as tangent vectors, such as v at point 
p G M in figure |3j Note here that whereas in the classical vector analysis the 
fields are mappings from the domain to a single vector space, the points of a 
manifold all have their distinct tangent spaces. 

The electric field is an object that yields the virtual electromotive force dU 
(virtual work up to charge) for every virtual translation v. In other words, the 
electric field is a functional field over the tangent spaces. Correspondingly, 
the electric flux is an object that yields the virtual flux for every pair of 
virtual translations, which define a virtual surface. That is, the electric flux 
is also a functional field over the tangent spaces. The solutions of quasi-static 
BVPs are such functional fields. 

2.5 Inner Product 

In the classical view we use a vector field to express the electric field. This 
possibility is due to Riesz representation theorem [12] [13], that enables us to 
represent the electric field as a vector field, once an inner product is avail- 
able. In the classical view the inner product is defined in the codomain 
vector space, but in the modern view each tangent space of M has its own 
inner product. The pointwise representation of the electric field functional 
by a vector E requires an inner product into each of them. At each point, 
with the inner product denoted by (•,■), the virtual emf dU{\) corresponding 
to the virtual translation vector v can be represented as dU{\) = (E,v). 
This representation makes the inner product an inherent part of the classical 
view: the vector E depends on the chosen inner product whereas dU does 
not. When the tangent spaces of a manifold are equipped with an inner 
product which varies smoothly from point to point, the manifold is called a 
Riemannian manifold [8j[9j[10j. The inner products on the tangent spaces of 
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Figure 4: Correspondences of tangent vectors under different charts. A tan- 
gent vector V at p G M corresponding to vectors dr/ at f{p) and dr^ at 

9{p)- 

M are usually called metric tensors. They induce a metric on M and also 
bring about the notion of angles between the tangent vectors [8]. 

For practical computations with tangent vectors, such as evaluations of 
metric tensors, we need the counterparts of the tangent vectors and the metric 
tensors on M". For each tangent vector v at p, a unique counterpart vector 
(the so-called push forward [llj) drj exists in M", induced by the chart /, see 
figure |4j We can use these counterpart vectors dr f to define metric tensors 
on M"' such that the tensors uniquely correspond to the ones on the tangent 
spaces of M: For each metric tensor (•, •)m given for M, there is a unique 
equivalent metric tensor (-, ■)f for M", that satisfies 

(dr/,dr'/)/ = (v,v')m (1) 

at all points and for all pairs dr/, dr'/ of M" and all pairs of corresponding 
tangent vectors v, v'. Equivalent metric tensors on MJ^ and M make / an 
isometry. However, the metric tensor chosen for M" need not be equivalent 
with the one on M, i.e. the charts need not be isometrics. 

^The metric tensors on the tangent spaces of M are selected such that the induced met- 
ric on M is the same as the metric M aheady has: the metric given by the measurements 
with some rigid body. 
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3 Equivalent Descriptions of BVP 



The theory of electromagnetics is independent of the chosen metric tensor 
and the chart, as stated in the introduction. In fact, Maxwell's equations do 
not need a metric tensor at all. They are meaningful on a differentiable man- 
ifold (not necessarily a Riemannian manifold) [E] [IS] [16] : the virtual works 
related to the displacements of a point charge are independent of the metric 
tensors and chosen charts. Naturally, the energy stored in a system does 
not depend on the chosen charts either. However, the constitutive equations 
need to be described to define the energy and there are at least two ways 
to proceed: First, a metric tensor is selected in order to represent the fields 
as vector fields, and constitutive equations are defined couple the two vector 
fields. Second way to proceed is to couple the functional fields directly,in 
which case the constitutive equations couple different types of fields. Al- 
though a metric tensor is not needed to construct the constitutive equations 
in this case, constitutive equations induce metric tensors on the manifold 
[16] . Yet in neither strategy does the energy itself depend on the selected 
or the induced metric tensor. Thus an electromagnetic BVP is independent 
of particular chart or metric tensor, and it is uniquely specified with mani- 
fold M, Maxwell's equations, the constitutive equations, and the boundary 
values. 

However, the material parameters, or rather the numeric values of per- 
mittivity, permeability, and conductivity depend on the chosen metric. Fur- 
thermore, to describe an abstract manifold concretely, we need to do that by 
charts, even though the manifold itself is independent of particular charts. 
Therefore, to pose the BVP, we need charts and metric tensors. When a 
chart and a metric tensor field are chosen, the numerical values of the ma- 
terial parameters can be specified. Then the material parameters are real- 
or matrix- valued fields defined on the chart. Thus we have a triplet {chart, 
metric tensor, material parameters} that describes a geometry of the BVP 
domain and the constitutive equation. Because we can choose any chart from 
the differentiable structure, and define any metric tensor on the chart, there 
are infinitely many geometries that we can use to pose the BVP. In fact, 
there is an equivalence class of these triplets {chart, metric tensor, material 
parameters}, that describe the same unique BVP. We call the equivalence 
relation producing these classes the material equivalence. 

We can now express a BVP on any chart: we need a triplet {chart, 
metric tensor, material parameters}. Maxwell's equations, and the boundary 
conditions. To express the BVP on any other chart, we only need to select 
another triplet, because Maxwell's equations are independent of the metric 
and of particular chart, as are the boundary conditions. Notice, that we 
can also keep the same chart and change the metric tensor and the material 



9 



parameters. These descriptions of the BVP have equivalent solutions; that 
is, in any description, the electric vector field corresponds to the same virtual 
emf functional field. We obtain a new view of EVPs from the material 
equivalence: a BVP can be seen as a material equivalence class together 
with Maxwell's equations and boundary conditions. This view suggests an 
approach to construct solver software systems: 

Solver for EVPs should allow the user to select any triplet from 
the material equivalence class. 

This way the user can use the triplet that is the most suitable for numerical 
solution. 

The benefit of choosing any triplet is that we can solve numerically a larger 
class of problems. Let us derive the material parameters for given chart and 
metric tensor, such that the material equivalence holds with another triplet 
that is known. 

4 Material Parameters for Given Chart And 
Metric Tensor 

To determine the material parameters for given chart and metric tensor, we 
require that the energy stored in a system and the virtual works related to 
the displacements of a point charge be invariant of our choice of chart and 
metric tensor. For simplicity, let us now focus on electrostatics, because other 
quasi-static cases can be treated conceptually in a similar fashion. 

Following the standard convention, we use dr ■ dr' to denote (dr, dr') in 
M". Also, for each metric tensor field • in M", there is a unique symmetric 
positive definite matrix field S such that dr • dr' = dr^S dr' holds pointwise 
for all dr and dr'. Let us assume we know a triplet {/,, e,}, consisting 
of chart fi of M, metric tensor -j on M" and the material parameters e^. 
Next, we choose chart fj from the difi^erentiable structure and use metric 
tensor -j on M" and then determine the material parameters ej for the triplet 
{fj, -j, Ej} such that the two triplets are equivalent. 

4.1 Invariance of Virtual Work 

The counterparts or the push forwards of a tangent vector v at p e M are 
drj and dr^ under charts fi and fj, respectively. The virtual emf dU{v) 
corresponding to the virtual translation v is Ej -j dr^, or Ej -j dr^, where -j 
and -j are metric tensors in R". These expressions must equal pointwise: 

E, dr,- = Ei -i dvi. (2) 
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Let fj o : ran(/i) C M" ^ M" be the transition map from the range of 
fi to the codomain of fj. Denoting by J the Jacobian matrix [17] of fj o /~\ 
we have 

drj = J dVi. (3) 

Substituting ^ for dr^ in (|2]), and relying on the invariance of the virtual 
work, we get 

EjS,Jdri = EjSidv, Vdr,,, (4) 

where Si and Sj are the matrix presentations of the metric tensors -j and -j. 
Because Q holds for every point and for every virtual translation dr^, we 
can express Ej in terms of E^ and matrices Si, Sj, and J as follows: 

Ei = S-^J^SjEj. (5) 

If we know either Ej or Ej, we know dU for every virtual translation, hence 
the solution of the BVP. 



4.2 Invariance of Energy 

The energy W stored in the electric field in fi{M) can be expressed as 

W = [ Ei-ieiEidvi= I EfSieiEidvi, (6) 

Jh(M) J MM) 

where the matrix contains the material parameters in the chart fi. Equa- 
tion induces an inner product 

{E,,Ei)= [ ElS,e,Eidv, (7) 

J MM) 

for the fields. Here E^ and E/ are any two vector fields. Even though there 
is a metric tensor in the expression of the inner product of the fields, the 
energy is independent of chosen metric tensor: the product Siei describes the 
constitutive equation, such that the value of energy is the same regardless 
of the choice of Si. This is what is meant by the claim that electromagnetic 
theory is independent of the chosen metric. Note that the inner product 
defined in ([T]) allows us to formulate the Galerkin method. 

For any electric field, the energy stored in the field should be invariant 
under the choice of chart, or 

/ EjSieiEidvi= [ EjSjejEjdvj (8) 

J MM) JfjiM) 
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must hold for any equivalent fields Ej and Ej. Using the change of variables 
theorem [17], we can write the left hand side of ([s]) in the codomain of fj as 



EjSi ei Ei dvi 



MM) 



Ef S'j ej Ej 



fj(M) 



J 



-1 



dvj, 



(9) 



where 



J- 



is the determinant of the inverse of the Jacobian matrix J. 



Substituting the expression given in equation (|5| for Ej in (joj), we get 

[ EjSjJS-^SieiS-\j^SjEj dvi = I EjSjejEjdvj. (10) 

JfAM) J fAM) 



Now, because this holds for any vector field Ej defined on fj{M), the matrix 
ej can be written in terms of matrices ej, Si, Sj, and J, as follows: 



€.j — J €{ J Sj 



(11) 



5 Equivalence of Numerical Solutions 

The solutions for equivalent triplets are equivalent from a mathematical point 
of view, and so are the numerical solutions: Next we define equivalent meshes 
and show that under equivalent meshes, the system matrices in FEM are 
identical up to the round-off errors of the floating point representation. 

First, assume a mesh in M. The mesh points in M have correspond- 
ing points in every codomain of the chart of every triplet; consequently, we 
get an equivalent mesh for every chart (an example of equivalent meshes in 
figure [5]) [^[^ Assuming real number arithmetics, equation ([s]) implies that 
the system matrices related to them are identical by construction. Thus the 
FEM solutions are also equivalent. However, in practice the system matrices 
are assembled and the numerical solutions are calculated with floating point 
numbers, and this leads to round-off errors. Thus the system matrices are 
identical only up to the round-off errors. 

Although the solutions are equivalent, in practice every EVP has triplets — 
sometimes including the standard parameterizations — that do not allow mesh 
generation, and thus no numerical solution. In these cases, mesh generation 
may well be possible with some other triplets. However, if the mesh quality 
criteria are not changed according to the chosen triplet, then there is usually 

^Notice that if the edges of mesh elements are straight in some chart, in an equivalent 
mesh of some other chart the edges can be curved. 

■^Elements that are long in one direction only can give good solutions if the field has 
little variations in that direction. 
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Figure 5: Example of equivalent meshes. Two equivalent meshes for the 
codomains of the charts / and g of figure |2} The mesh in the codomain 
of the standard parameterization / has long elements in the middle, which 
may imply a poor quality. However, it is not always possible to generate any 
reasonable mesh for the standard parameterizations. 

a price to pay: the mesh quality may be poor, which causes numerical error 
in the solution. For example some elements may be very long in the standard 
parameterization as seen in top of figure |5| 

6 Applications 

We now have all the necessary tools to express a EVP, attempt a solution, 
and interpret it with any triplet, and next suggest how to take advantage 
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of this possibility in practice. First, we show how a chart can be changed 
when the metric tensor in M" is fixed, which is usually the case. Second, 
we suggest a scheme to use multiple charts in mesh generation and problem 
solution. Third application is open boundary problems which are solved by 
selecting a triplet in which the problem domain in the chart is bounded. 
Finally, we discuss how to model motion with only one mesh by changing the 
metric. 



6.1 Reparameterization Under Fixed Chart Metric 

Because the numerical values of the material parameters are always given 
in literature in some length unit system, it is easiest to describe the con- 
stitutive relations with isometric standard parameterizations. Consequently, 
the programmer of a finite element software system often assumes that the 
user applies only isometric standard parameterizations; that is, the user can 
control only the unit of length, and rotate and translate the model. Further- 
more, due to the isometry assumption, the metric tensor is implicitly fixed 
such that it induces the Euclidean 2-norm. Therefore the matrix representing 
the metric tensor in the intrinsic coordinate system of the solver is implicitly 
hardwired to be the identity matrix. 

Let us assume that we know how to express a BVP in some isometric 
standard parameterization /, that is, we know the triplet {/, Sf, e/}. How- 
ever, instead of using / directly, we would like to reparameterize with chart 
g, which is not a standard parameterization, to avoid fioating point problems 
in numerical modeling. Because of the same hardwired metric tensor in the 
codomains of / and g, the distances between any two coordinates (tuples of 
real numbers) do not change when changing chart from / to g. However, 
a change of chart makes the coordinates correspond with different points in 
the manifold M, and thus the change of chart in this case implies a change 
of metric, hence the chart g is not isometry. 

To set up our BVP correctly under g, we must impose the material equiv- 
alence on g. For this we can use only the material parameters eg, because 
g is chosen by the modeller, and the metric tensor is hardwired by the pro- 



grammer. Equation (11) shows that the material parameters eg require in- 
formation on g relative to /. Yet everything related to /, including / it- 
self, exists only in the modeller's mind. We are solving the BVP using a 
triplet {g, Sg, eg}, and we interpret this description of the BVP under g as 
being equivalent to that under /. However, the software programmer who 
assumes only standard parameterizations and has no information about /, 
might think of (7 as a standard parameterization and very likely interpret the 
problem geometry very differently. In this sense, we would be "misusing" the 
software. 
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To impose the material equivalence, we need the Jacobian between the 
charts and material parameters e/. After substituting I for Sf and Sg in 
(11), we get the following equivalent material parameters: 



J 



-JefJ^ 



(12) 



where the identity 



J 



J 



-1 



1 is used (Henrotte et al. yjj arrived at a 



similar expression but one of their motivations was open boundary problems, 



see section 6.3). It is easy to see from (12) that the material parameters in g 



may not be a scalar multiple of the identity matrix, even if they were in the 
chart /, a point expressed in [16] from the opposite point of view. 

The solution of the problem are fields and Dp. However, to present the 
fields in the standard parameterization /, where they may be more intuitive, 
we must use ([S]) and = e/Ej to obtain the equivalent fields. 

Remark 1 cannot he transformed into D / with the same transformation 
as Eg is transformed into Ef in That is, the electric field intensity E 
and the electric flux density D transform differently between charts. To get 
the Df, one can use ^ to transform Eg to Ej and then use the constitutive 
law Df = 6/ E J . 

The above method is fully realizable with any mesh generator and solver 
software that allows for pointwise definition of material parameters as a ma- 
trix. However, if the solver does not allow for pointwise definition, only charts 
g with piecewise affine transition maps from the standard parameterization 
/ can be used. 



6.2 Modelling With Multiple Charts 

Previously we had two charts that both cover the whole domain. However, 
manifold can be covered with multiple charts, such that only a part of the 
manifold is covered by each chart. All that is required is that every point 
of the manifold belongs to the domain of some chart and that all the charts 
belong to the same differentiable structure. The use of multiple charts can 
benefit numerical analysis, for example in the case of figure |6] First of all, the 
mesh generation gets easier, because the detailed regions could be covered 
separately with their own charts, with the origins located such that there are 
as many floating point numbers available as possible. Second, for the same 
reason, the calculations of the inner products for the system matrix are more 
accurate. 
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Figure 6: Multiple charts. Top: Standard parameterization, with three rect- 
angular regions. Bottom: Three charts, that each cover one of the rectangular 
regions, but the origins are moved and scales are changed. 



Even though multiple charts are used, we still assume here that the mani- 
fold can be covered with one M"" chart. This makes it possible to partition the 
manifold into separate domains such that the domains only overlap on their 
boundary. Furthermore, this guarantees that the manifold can be embedded 
into the ra-dimensional Euclidean space. Next, based on these assumptions, 
we suggest an implementation of multiple charts into software systems. 

First, the user implicitly decides on one isometric standard parameteri- 
zation, which covers the whole domain. Let us call it the universal chart, 
because it works as the reference to the other charts, and the material pa- 
rameters are given on it. However, the universal chart is not used in the 
calculation or in the mesh generation, but these are done in separate charts 
that cover only parts of the manifold: the user partitions the universal chart 
into multiple regions, such that the regions only overlap on their boundary 
(see figure [gJ. Then the user gives new charts for these regions and gives the 
transition maps to the universal chart. In practice the user only gives the 
separate charts and the transition maps to the software system, and then the 
universal chart is constructed from these by the software. It is also possi- 
ble that the universal chart is not constructed at all in the software system 
but exists only implicitly. The metric tensors for these separate charts could 
be fixed or given by the user. The material parameters for the charts are 
given as if they are standard parameterizations, and the software calculates 
automatically the material parameters to be used in the calculations. This 
is possible because the metric tensors are known, the transition maps, 

whose Jacobians can be automatically calculated. 

Finally, in order to attain compatibility of meshes in different charts, the 
meshes have to agree at the boundaries of the regions. The mesh generation 
could proceed as follows: First a mesh is generated on one of the regions with 
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Figure 7: Open boundary problem. Top: The problem domain (manifold) 
consist of an device over the ground. Bottom: The codomains of two charts, 
/ and g. The regions inside of the rectangles denoted with dashed line cor- 
respond to the same part of the domain, but in g the rest of the space is 
between the rectangles. 

the chart that covers it. Then the coordinates of the nodes on the common 
boundary with a second region are mapped with the transition map to the 
chart that covers the second region, and the rest of the mesh in the second 
region can be generated. The assembly of the system matrix can be done as 
usually. 

6.3 Open Boundary Problems 

The open boundary problems are often solved by truncating the domain and 
setting the fields zero at an artificial distant boundary. This usually gives 
a good enough approximation for the fields near the sources, because the 
fields tend to zero quite fast as the distance from the source increases. The 
distance of the artificial boundary from the sources and regions of importance 
is decided by the modeller. However, the error from this truncation is hard 
to estimate and the number of elements needed to cover all that empty space 
can be large. 

To overcome these problems, many methods have been proposed over the 
years: The so called "ballooning method" [19], in which the true distance 
of the boundary is put very far away with a thin layer of special elements. 
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Another method uses infinite elements [20], which are special decaying basis 
functions used in FEM. Yet another method couples the FEM with analytical 
solutions as in [21]. Finally, the transformation (or shell transformation) 
methods presented for example in [18] [22] put the infinite boundary to a finite 
distance with the help of a suitable change of coordinates. However, many of 
these methods are hard to implement into most production solvers, because 
they require modification of the solver codes: For example, the modeller 
needs to define new basis functions or needs to give pointwise Jacobians of 
the change of coordinates. 

The triplet approach of EVPs suggest a strategy to solve open boundary 
problems: Just select a suitable triplet, such that in the chart the infinity 
boundary is at a finite distance, such as the chart g in the figure [7} This 
is effectively identical to many of the above methods. For example, in the 
transformations method the Jacobian of the change of coordinates is imple- 
mented in the software code whereas in the triplet approach the Jacobian is 
included in the material parameters given by the user. 

The triplet should be selected such that a small region including the 
source regions and other regions of interest is modelled as usually, but the 
complement of this region is mapped to a small bounded outer region as in 
the chart g in the figure [7| The elements nearest the infinity boundary of the 
domain are such that the points at the boundary are mapped to the infinity by 
the transition map form g to the corresponding standard parameterization. 
A suitable transition mapping could be the following: The region outside the 
unit disc is mapped to the ring with exterior radius two, such that the ring 
surrounds the unit disc and the boundary of the disc is the inner boundary 
of the ring, see figure [8] A point outside the unit disc at the distance r from 
the center of the unit disc is mapped to a point at the distance R from the 
center, such that the distances have the following relation: 

R = 2--. 

r 

This method can be used with almost any solver: Only the material pa- 
rameters have to be changed accordingly, which is a responsibility of the 
modeller. However, the procedure is easy to automate: the user only spec- 
ifies the regions where the outer boundary is mapped to the infinity by the 
transition map from g to the standard parameterization, and the software 
system would take care of the rest. The system could have built-in such a 
mapping or mappings which adapt according to the coordinates of the in- 
finity boundary and the boundary of the infinite regions. If necessary, the 
user gives the correspondence between these boundaries and then the sys- 
tem modifies the user-given material parameters according to the adapted 
mapping. 
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Figure 8: Mapping of infinite domain to finite domain. The region outside 
of the unit disc in the left hand side of the figure is mapped to ring of radius 
two in the right hand side of the figure. A point at distance r from the center 
of unit disc is mapped to a point at distance R from the center. 

6.4 Modelling Motion With Single Chart 

When modelhng motion, we have to solve a sequence of BVPs whose only 
difference is that some object changes its position relative to the other ob- 
jects. For an example, consider the force calculations of an electromagnet, 
see figure |9| The magnet attracts the load, which then moves upwards. One 
solution is to make whole model for each position of the moving object. For 
each position, we have to define a new geometry and generate a new mesh 
for it, and furthermore assembly the whole system matrix again. This is very 
time-consuming. 

The triplet approach can speed up the process: Keep the chart and the 
material parameters fixed and change the metric tensor. This changes the 
geometry and takes into account the motion: When a rigid object moves in a 
fluid (air) such that the topology does not change, i.e. the moving object does 
not touch any other material than the fluid, then the underlying differentiable 
manifold of the domain does not change (this underlying manifold is often 
referred to as material manifold). However, the Riemannian manifold will 
change, because the metric tensor changes. This means that the distances 
between points of the manifold are changed. Thus we may select a single 
chart g from the differentiable structure, to cover all the different Riemannian 
manifolds. 

The metric tensor for the selected chart g again depends on the material 
equivalence. This requires that some standard parameterization / for each 
Riemannian manifold is known, when the corresponding Jacobians are also 
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Figure 9: Modelling motion. Electromagnet: coil with current /, ferromag- 
netic core and load, which is attracted upwards. 



known. Then the matrix Sg representing the metric tensor on g can be solved 



from equation (11): Assume the material parameter e/ to be a scalar and 
require the material parameters to be equal, that is eg = e/. Then Sg in 
terms of the Jacobian can be solved: 



Sg 



J 



j-^.r\ (13) 



Now, for the example in the case of figure [9} the metric tensor only changes 



in the air between the magnet and the load according to the equation (13). 
Because the mesh remains the same, only the elements in the system matrix 
corresponding to the mesh elements in this air region will change. Thus 
only a partial re-assembly of system matrix is needed for each new step. 
Furthermore, in the case of iterative solvers, a preconditioner used for one 
step can be used effectively for many successive steps, because the changes 
in the system matrix can be very small. Also the previous solution may be a 
good initial guess. Thus, all these things — only one mesh, partial re-assembly 
of system matrix, the preconditioner, and the initial guess — can speed up the 
solution process substantially. 

Many solvers do not offer the possibility to change the metric tensor, 
yet the modelling of motion with a single mesh is still possible: we can 
equivalently change the material parameters, see section 6.1. One needs to 
know a standard parameterization for each Riemannian manifold and then 



use the material parameters given by the equation (12). Finally, all things 
said about modelling motion apply also for modelling the geometric effects 
of deformations, such as the effects of magnetostriction. 
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Figure 10: Power line in a standard parameterization. 




Figure 11: The chart used in the power hne problem. 

7 Numerical Example 

We demonstrate a combination of proposed applications: an open bound- 
ary problem (section 6.3) with reparameterization under fixed chart metric 
(section 6.1). The problem to solve is the electric potential of a three-phase 
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Figure 12: Isovalue lines of potential of the power line shown in the chart 
used in the calculations. 



high-voltage line, a 3D-Laplace's problem. Figure 10 shows a standard pa- 
rameterization of the domain. The domain is the half-space above the ground 
apart from the line itself. The infinity boundary is mapped to a finite dis- 
tance from the line, as explained in section 6.3. There is also great variation 
in the scale of the details of the domain: the length of the lines is order 
of hundreds of meters, whereas the smallest dimensions are order of cen- 
timeters. Based on experience, we can say that the mesh generation with 
the standard parameterization is difficult or even impossible. To avoid mesh 
generation problems, the lines are scaled down, as is the height of the pillars. 
The model used in the calculations is shown in figure 11 The result is shown 
in figure 12 The figures and the calculations are produced with GetDP [23] 
and Gmsh 



8 Conclusion 

Most solver software systems for quasi-static electromagnetic boundary value 
problems are based on the formulation of electromagnetics with classical vec- 
tor analysis. This has many restrictions for numerical solution, because many 
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mathematical structures have hmited freedom or they are completely fixed. 
For example, the fields are presented as mappings from some global coor- 
dinate system to a single vector space and the metrics for the coordinate 
systems and the vector space are fixed. The restrictions limit the possible 
coordinate systems that can be used. 

The formulation of electromagnetics with manifolds and differential ge- 
ometry is less restrictive and thus can help in a numerical solution. This 
view exposes the independence of the electromagnetic theory on the chosen 
coordinates or metric. Based on this it is shown that there is an equivalence 
class of triplets {coordinate system, metric, material parameters} for each 
boundary value problem. A triplet describes a geometry of the problem do- 
main and the constitutive equations which are needed to pose the problem in 
the coordinate system. We propose that the solvers should allow to use any 
triplet, that is, the modeller should be able to choose any coordinate system 
and metric at will. These possibilities come built-in in differential geometry, 
whereas they are not self-evident in classical vector analysis. We discuss how 
this can be done. It is not hard to apply some of the triplet strategies in 
most production solvers, because not all of them require modifications to 
the code. We show that despite the freedom in the choice of the triplet, 
the FEM system matrices of a given BVP, assembled for a mesh represented 
under different triplets, arc identical. 

We show how the triplet approach can be exploited in applications: Po- 
tentially prohibitive difficulties of mesh generation in problem domains that 
have large scale variation can be alleviated by a suitable choice of coordinate 
system. This extends the class of solvable problems, although possibly at the 
expense of mesh quality. Open boundary problems can be solved by choosing 
a convenient triplet for the solution. Solution of problems involving motion 
and deformations can be accelerated when a single mesh is used on every 
time step and only the metric is altered. 
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